library(SingleCellExperiment)
library(BiocParallel)
library(geneBasisR)
library(splatter)
library(scran)
library(miloR)
library(miloDE)
library(stringr)
library(ggplot2)
library(ggpubr)
library(viridis)
library(MetBrewer)
library(dplyr)
library(tibble)
ncores = 24
mcparam = MulticoreParam(workers = ncores)
register(mcparam)
set.seed(32)
root.dir = "/nfs/research/marioni/alsu/hubmap_metaRef/"
source(paste0(root.dir , "miloDE_analysis/core_functions.R"))
figures.dir = paste0(root.dir , "figures/da_de/simulations_splatter/random/")
# generated simulation data
sces = readRDS(file = paste0(root.dir , "data/processed/da_de/simulations_random/sims.Rds"))
simulated_genes = readRDS(file = paste0(root.dir , "data/processed/da_de/simulations_random/simulated_genes.Rds"))
stat_da = readRDS(file = paste0(root.dir , "data/processed/da_de/simulations_random/stat_da.Rds"))
stat_de = readRDS(file = paste0(root.dir , "data/processed/da_de/simulations_random/stat_de.Rds"))
stat_combined = readRDS(file = paste0(root.dir , "data/processed/da_de/simulations_random/stat_combined_all_genes.Rds"))
stat_combined_reduced = unique(stat_combined[, c("id_sim" , "id_nhood_assignment" , "reducedDim_name", "frac_hoods_DA_001", "frac_hoods_DA_005", "frac_hoods_DA_01" , "n_markers" , "n_not_markers" , "round_sim" , "order" , "k" , "seed")] )
stat_combined_reduced$order_k = paste0(stat_combined_reduced$order , "_" , stat_combined_reduced$k)
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.hvgs"] = "Supervised"
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.all"] = "Unsupervised"
get_plot = function(var){
p = ggplot(stat_combined_reduced , aes(x = factor(n_markers) , y = factor(n_not_markers) , fill = stat_combined_reduced[,var])) +
geom_tile() +
facet_wrap(~reducedDim_name + round_sim + order_k , ncol = 12) +
scale_fill_viridis(discrete = F , name = "Fraction of DA nhoods") +
theme_bw() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) +
labs(x = "# of altered genes contributing to joint variance" , y = "# of altered genes contributing to case variance") +
theme(legend.position = "top")
p
return(p)
}
p = get_plot("frac_hoods_DA_01")
ggsave(filename = paste0(figures.dir, "da/split_all__fdr_01", ".png"), plot = p, width = 20, height = 6)
p = get_plot("frac_hoods_DA_005")
ggsave(filename = paste0(figures.dir, "da/split_all__fdr_005", ".png"), plot = p, width = 20, height = 6)
p = get_plot("frac_hoods_DA_001")
ggsave(filename = paste0(figures.dir, "da/split_all__fdr_001", ".png"), plot = p, width = 20, height = 6)
stat_combined_reduced = unique(stat_combined[, c("id_sim" , "id_nhood_assignment" , "reducedDim_name", "frac_hoods_DA_001", "frac_hoods_DA_005", "frac_hoods_DA_01" , "n_markers" , "n_not_markers" , "round_sim" , "order" , "k" , "seed")] )
stat_combined_reduced$order_k = paste0(stat_combined_reduced$order , "_" , stat_combined_reduced$k)
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.hvgs"] = "Supervised"
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.all"] = "Unsupervised"
stat_combined_reduced = as.data.frame(stat_combined_reduced %>% group_by(id_sim , reducedDim_name, frac_hoods_DA_001, frac_hoods_DA_005, frac_hoods_DA_01 , n_markers , n_not_markers , round_sim) %>% dplyr::summarise(frac_hoods_DA_001 = mean(frac_hoods_DA_001 , na.rm = T) , frac_hoods_DA_005 = mean(frac_hoods_DA_005 , na.rm = T) , frac_hoods_DA_01 = mean(frac_hoods_DA_01 , na.rm = T)) )
get_plot = function(var){
p = ggplot(stat_combined_reduced , aes(x = factor(n_markers) , y = factor(n_not_markers) , fill = stat_combined_reduced[,var])) +
geom_tile() +
facet_wrap(~reducedDim_name + round_sim) +
scale_fill_viridis(discrete = F , name = "Fraction of DA nhoods") +
theme_bw() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) +
labs(x = "# of altered genes contributing to joint variance" , y = "# of altered genes contributing to case variance") +
theme(legend.position = "top")
p
return(p)
}
p = get_plot("frac_hoods_DA_01")
ggsave(filename = paste0(figures.dir, "da/groupped_sim__fdr_01", ".png"), plot = p, width = 8, height = 6)
p = get_plot("frac_hoods_DA_005")
ggsave(filename = paste0(figures.dir, "da/groupped_sim__fdr_005", ".png"), plot = p, width = 8, height = 6)
p = get_plot("frac_hoods_DA_001")
ggsave(filename = paste0(figures.dir, "da/groupped_sim__fdr_001", ".png"), plot = p, width = 8, height = 6)
stat_combined_reduced = unique(stat_combined[, c("id_sim" , "id_nhood_assignment" , "reducedDim_name", "frac_hoods_DA_001", "frac_hoods_DA_005", "frac_hoods_DA_01" , "n_markers" , "n_not_markers" , "round_sim" , "order" , "k" , "seed")] )
stat_combined_reduced$order_k = paste0(stat_combined_reduced$order , "_" , stat_combined_reduced$k)
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.hvgs"] = "Supervised"
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.all"] = "Unsupervised"
stat_combined_reduced = as.data.frame(stat_combined_reduced %>% group_by(reducedDim_name, frac_hoods_DA_001, frac_hoods_DA_005, frac_hoods_DA_01 , n_markers , n_not_markers ) %>% dplyr::summarise(frac_hoods_DA_001 = mean(frac_hoods_DA_001 , na.rm = T) , frac_hoods_DA_005 = mean(frac_hoods_DA_005 , na.rm = T) , frac_hoods_DA_01 = mean(frac_hoods_DA_01 , na.rm = T)) )
get_plot = function(var){
p = ggplot(stat_combined_reduced , aes(x = factor(n_markers) , y = factor(n_not_markers) , fill = stat_combined_reduced[,var])) +
geom_tile() +
facet_wrap(~reducedDim_name) +
scale_fill_viridis(discrete = F , name = "Fraction of DA nhoods") +
theme_bw() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) +
labs(x = "# of altered genes contributing to joint variance" , y = "# of altered genes\ncontributing to case variance") +
theme(legend.position = "top")
p
return(p)
}
p = get_plot("frac_hoods_DA_01")
ggsave(filename = paste0(figures.dir, "da/groupped_all__fdr_01", ".png"), plot = p, width = 5, height = 3)
p = get_plot("frac_hoods_DA_005")
ggsave(filename = paste0(figures.dir, "da/groupped_all__fdr_005", ".png"), plot = p, width = 5, height = 3)
p = get_plot("frac_hoods_DA_001")
ggsave(filename = paste0(figures.dir, "da/groupped_all__fdr_001", ".png"), plot = p, width = 5, height = 3)
stat_combined_reduced = unique(stat_combined[, c("id_sim" , "id_nhood_assignment" , "reducedDim_name", "gene_type", "frac_hoods_DE_001", "frac_hoods_DE_005", "frac_hoods_DE_01" , "n_markers" , "n_not_markers" , "round_sim" , "order" , "k" , "seed")] )
stat_combined_reduced$order_k = paste0(stat_combined_reduced$order , "_" , stat_combined_reduced$k)
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.hvgs"] = "Supervised"
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.all"] = "Unsupervised"
get_plot = function(var , gene_type){
current.stat = stat_combined_reduced[stat_combined_reduced$gene_type == gene_type , ]
p = ggplot(current.stat , aes(x = factor(n_markers) , y = factor(n_not_markers) , fill = current.stat[,var])) +
geom_tile() +
facet_wrap(~reducedDim_name + round_sim + order_k , ncol = 12) +
scale_fill_viridis(discrete = F , name = "Fraction of DE nhoods") +
theme_bw() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) +
labs(x = "# of altered genes contributing to joint variance" , y = "# of altered genes contributing to case variance") +
theme(legend.position = "top")
p
return(p)
}
anno = expand.grid(gene_type = c("marker" , "not_marker") , var = c("01" , "005" , "001"))
plots = lapply(1:nrow(anno) , function(i){
p = get_plot(paste0("frac_hoods_DE_" ,as.character(anno$var[i])) , as.character(anno$gene_type[i]))
ggsave(filename = paste0(figures.dir, "de/" , as.character(anno$gene_type[i]), "/split_all__fdr_" , as.character(anno$var[i]), ".png"), plot = p, width = 20, height = 6)
})
stat_combined_reduced = unique(stat_combined[, c("id_sim" , "id_nhood_assignment" , "reducedDim_name", "gene_type", "frac_hoods_DE_001", "frac_hoods_DE_005", "frac_hoods_DE_01" , "n_markers" , "n_not_markers" , "round_sim" , "order" , "k" , "seed")] )
stat_combined_reduced$order_k = paste0(stat_combined_reduced$order , "_" , stat_combined_reduced$k)
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.hvgs"] = "Supervised"
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.all"] = "Unsupervised"
stat_combined_reduced = as.data.frame(stat_combined_reduced %>% group_by(id_sim , reducedDim_name, gene_type , n_markers , n_not_markers , round_sim) %>% dplyr::summarise(frac_hoods_DE_001 = mean(frac_hoods_DE_001 , na.rm = T) , frac_hoods_DE_005 = mean(frac_hoods_DE_005 , na.rm = T) , frac_hoods_DE_01 = mean(frac_hoods_DE_01 , na.rm = T)) )
get_plot = function(var , gene_type){
current.stat = stat_combined_reduced[stat_combined_reduced$gene_type == gene_type , ]
p = ggplot(current.stat , aes(x = factor(n_markers) , y = factor(n_not_markers) , fill = current.stat[,var])) +
geom_tile() +
facet_wrap(~reducedDim_name + round_sim ) +
scale_fill_viridis(discrete = F , name = "Fraction of DE nhoods") +
theme_bw() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) +
labs(x = "# of altered genes contributing to joint variance" , y = "# of altered genes contributing to case variance") +
theme(legend.position = "top")
p
return(p)
}
anno = expand.grid(gene_type = c("marker" , "not_marker") , var = c("01" , "005" , "001"))
plots = lapply(1:nrow(anno) , function(i){
p = get_plot(paste0("frac_hoods_DE_" ,as.character(anno$var[i])) , as.character(anno$gene_type[i]))
ggsave(filename = paste0(figures.dir, "de/" , as.character(anno$gene_type[i]), "/groupped_sim__fdr_" , as.character(anno$var[i]), ".png"), plot = p, width = 8, height = 6)
})
stat_combined_reduced = unique(stat_combined[, c("id_sim" , "id_nhood_assignment" , "reducedDim_name", "gene_type", "frac_hoods_DE_001", "frac_hoods_DE_005", "frac_hoods_DE_01" , "n_markers" , "n_not_markers" , "round_sim" , "order" , "k" , "seed")] )
stat_combined_reduced$order_k = paste0(stat_combined_reduced$order , "_" , stat_combined_reduced$k)
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.hvgs"] = "Supervised"
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.all"] = "Unsupervised"
stat_combined_reduced = as.data.frame(stat_combined_reduced %>% group_by(reducedDim_name, gene_type , n_markers , n_not_markers , round_sim) %>% dplyr::summarise(frac_hoods_DE_001 = mean(frac_hoods_DE_001 , na.rm = T) , frac_hoods_DE_005 = mean(frac_hoods_DE_005 , na.rm = T) , frac_hoods_DE_01 = mean(frac_hoods_DE_01 , na.rm = T)) )
get_plot = function(var , gene_type){
current.stat = stat_combined_reduced[stat_combined_reduced$gene_type == gene_type , ]
p = ggplot(current.stat , aes(x = factor(n_markers) , y = factor(n_not_markers) , fill = current.stat[,var])) +
geom_tile() +
facet_wrap(~reducedDim_name ) +
scale_fill_viridis(discrete = F , name = "Fraction of DE nhoods" ) +
theme_bw() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) +
labs(x = "# of altered genes contributing to joint variance" , y = "# of altered genes\ncontributing to case variance") +
theme(legend.position = "top")
p
return(p)
}
anno = expand.grid(gene_type = c("marker" , "not_marker") , var = c("01" , "005" , "001"))
plots = lapply(1:nrow(anno) , function(i){
p = get_plot(paste0("frac_hoods_DE_" ,as.character(anno$var[i])) , as.character(anno$gene_type[i]))
ggsave(filename = paste0(figures.dir, "de/" , as.character(anno$gene_type[i]), "/groupped_all__fdr_" , as.character(anno$var[i]), ".png"), plot = p, width = 5, height = 3)
})
stat_combined_reduced = unique(stat_combined[, c("id_sim" , "id_nhood_assignment" , "reducedDim_name", "gene_type", "frac_hoods_DE_not_DA_001", "frac_hoods_DE_not_DA_005", "frac_hoods_DE_not_DA_01" , "n_markers" , "n_not_markers" , "round_sim" , "order" , "k" , "seed")] )
stat_combined_reduced$order_k = paste0(stat_combined_reduced$order , "_" , stat_combined_reduced$k)
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.hvgs"] = "Supervised"
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.all"] = "Unsupervised"
get_plot = function(var , gene_type){
current.stat = stat_combined_reduced[stat_combined_reduced$gene_type == gene_type , ]
p = ggplot(current.stat , aes(x = factor(n_markers) , y = factor(n_not_markers) , fill = current.stat[,var])) +
geom_tile() +
facet_wrap(~reducedDim_name + round_sim + order_k , ncol = 12) +
scale_fill_viridis(discrete = F , name = "Fraction of DE nhoods\n(in not DA nhoods)") +
theme_bw() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) +
labs(x = "# of altered genes contributing to joint variance" , y = "# of altered genes contributing to case variance") +
theme(legend.position = "top")
p
return(p)
}
anno = expand.grid(gene_type = c("marker" , "not_marker") , var = c("01" , "005" , "001"))
plots = lapply(1:nrow(anno) , function(i){
p = get_plot(paste0("frac_hoods_DE_not_DA_" ,as.character(anno$var[i])) , as.character(anno$gene_type[i]))
ggsave(filename = paste0(figures.dir, "de_in_not_da/" , as.character(anno$gene_type[i]), "/split_all__fdr_" , as.character(anno$var[i]), ".png"), plot = p, width = 20, height = 6)
})
stat_combined_reduced = unique(stat_combined[, c("id_sim" , "id_nhood_assignment" , "reducedDim_name", "gene_type", "frac_hoods_DE_not_DA_001", "frac_hoods_DE_not_DA_005", "frac_hoods_DE_not_DA_01" , "n_markers" , "n_not_markers" , "round_sim" , "order" , "k" , "seed")] )
stat_combined_reduced$order_k = paste0(stat_combined_reduced$order , "_" , stat_combined_reduced$k)
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.hvgs"] = "Supervised"
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.all"] = "Unsupervised"
stat_combined_reduced = as.data.frame(stat_combined_reduced %>% group_by(id_sim , reducedDim_name, gene_type , n_markers , n_not_markers , round_sim) %>% dplyr::summarise(frac_hoods_DE_001 = mean(frac_hoods_DE_not_DA_001 , na.rm = T) , frac_hoods_DE_005 = mean(frac_hoods_DE_not_DA_005 , na.rm = T) , frac_hoods_DE_01 = mean(frac_hoods_DE_not_DA_01 , na.rm = T)) )
get_plot = function(var , gene_type){
current.stat = stat_combined_reduced[stat_combined_reduced$gene_type == gene_type , ]
p = ggplot(current.stat , aes(x = factor(n_markers) , y = factor(n_not_markers) , fill = current.stat[,var])) +
geom_tile() +
facet_wrap(~reducedDim_name + round_sim ) +
scale_fill_viridis(discrete = F , name = "Fraction of DE nhoods\n(in not DA nhoods)") +
theme_bw() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) +
labs(x = "# of altered genes contributing to joint variance" , y = "# of altered genes contributing to case variance") +
theme(legend.position = "top")
p
return(p)
}
anno = expand.grid(gene_type = c("marker" , "not_marker") , var = c("01" , "005" , "001"))
plots = lapply(1:nrow(anno) , function(i){
p = get_plot(paste0("frac_hoods_DE_" ,as.character(anno$var[i])) , as.character(anno$gene_type[i]))
ggsave(filename = paste0(figures.dir, "de_in_not_da/" , as.character(anno$gene_type[i]), "/groupped_sim__fdr_" , as.character(anno$var[i]), ".png"), plot = p, width = 8, height = 6)
})
stat_combined_reduced = unique(stat_combined[, c("id_sim" , "id_nhood_assignment" , "reducedDim_name", "gene_type", "frac_hoods_DE_not_DA_001", "frac_hoods_DE_not_DA_005", "frac_hoods_DE_not_DA_01" , "n_markers" , "n_not_markers" , "round_sim" , "order" , "k" , "seed")] )
stat_combined_reduced$order_k = paste0(stat_combined_reduced$order , "_" , stat_combined_reduced$k)
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.hvgs"] = "Supervised"
stat_combined_reduced$reducedDim_name[stat_combined_reduced$reducedDim_name == "pca.all"] = "Unsupervised"
stat_combined_reduced = as.data.frame(stat_combined_reduced %>% group_by(reducedDim_name, gene_type , n_markers , n_not_markers , round_sim) %>% dplyr::summarise(frac_hoods_DE_001 = mean(frac_hoods_DE_not_DA_001 , na.rm = T) , frac_hoods_DE_005 = mean(frac_hoods_DE_not_DA_005 , na.rm = T) , frac_hoods_DE_01 = mean(frac_hoods_DE_not_DA_01 , na.rm = T)) )
get_plot = function(var , gene_type){
current.stat = stat_combined_reduced[stat_combined_reduced$gene_type == gene_type , ]
p = ggplot(current.stat , aes(x = factor(n_markers) , y = factor(n_not_markers) , fill = current.stat[,var])) +
geom_tile() +
facet_wrap(~reducedDim_name ) +
scale_fill_viridis(discrete = F , name = "Fraction of DE nhoods" ) +
theme_bw() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) +
labs(x = "# of altered genes contributing to joint variance" , y = "# of altered genes\ncontributing to case variance") +
theme(legend.position = "top")
p
return(p)
}
anno = expand.grid(gene_type = c("marker" , "not_marker") , var = c("01" , "005" , "001"))
plots = lapply(1:nrow(anno) , function(i){
p = get_plot(paste0("frac_hoods_DE_" ,as.character(anno$var[i])) , as.character(anno$gene_type[i]))
ggsave(filename = paste0(figures.dir, "de_in_not_da/" , as.character(anno$gene_type[i]), "/groupped_all__fdr_" , as.character(anno$var[i]), ".png"), plot = p, width = 5, height = 3)
})
@@ Mb later, lets see if we need it
cols = met.brewer("Egypt",4)
cols = cols[c(3:4)]
get_plot = function(id , red_dim){
sce = sces[[id]]
umaps = as.data.frame( scater::calculateUMAP(t(reducedDim(sce , red_dim))) )
reducedDim(sce , "umap") = umaps
meta = cbind(as.data.frame(colData(sce)) , as.data.frame(reducedDim(sce , "umap")))
meta$Group = as.character(meta$Group)
meta$Group[meta$Group == "Group1"] = "Cell type 1"
meta$Group[meta$Group == "Group2"] = "Cell type 2"
meta$type[meta$type == "cond_A"] = "Condition A"
meta$type[meta$type == "cond_B"] = "Condition B"
meta = meta[sample(nrow(meta)) , ]
#title = names(sces)[id]
#title = str_remove(title , "__n_not_markers_0")
p = ggplot(meta , aes(x = V1 , y = V2 , col = factor(Group))) +
geom_point(alpha = .5 , size = .5) +
theme_bw() +
facet_wrap(~type, nrow = 2) +
labs(x = "UMAP-1" , y = "UMAP-2") +
scale_color_manual(values = cols , name = "Cell type")
p
return(p)
}
id = 1
p = get_plot(id, "pca.hvgs")
p
ggsave(filename = paste0(figures.dir, "umaps/hvgs/" , names(sces)[id], ".png"), plot = p, width = 4, height = 5)
p = get_plot(id, "pca.all")
p
ggsave(filename = paste0(figures.dir, "umaps/all/" , names(sces)[id], ".png"), plot = p, width = 4, height = 5)
id = 25
p = get_plot(id, "pca.hvgs")
p
ggsave(filename = paste0(figures.dir, "umaps/hvgs/" , names(sces)[id], ".png"), plot = p, width = 4, height = 5)
p = get_plot(id, "pca.all")
p
ggsave(filename = paste0(figures.dir, "umaps/all/" , names(sces)[id], ".png"), plot = p, width = 4, height = 5)
id = 361
p = get_plot(id, "pca.hvgs")
p
ggsave(filename = paste0(figures.dir, "umaps/hvgs/" , names(sces)[id], ".png"), plot = p, width = 4, height = 5)
p = get_plot(id, "pca.all")
p
ggsave(filename = paste0(figures.dir, "umaps/all/" , names(sces)[id], ".png"), plot = p, width = 4, height = 5)
id = 385
p = get_plot(id, "pca.hvgs")
p
ggsave(filename = paste0(figures.dir, "umaps/hvgs/" , names(sces)[id], ".png"), plot = p, width = 4, height = 5)
p = get_plot(id, "pca.all")
p
ggsave(filename = paste0(figures.dir, "umaps/all/" , names(sces)[id], ".png"), plot = p, width = 4, height = 5)
get_all_plots = function(round_sim , red_dim){
ids = sapply(names(sces) , function(str) str_detect(str , paste0(round_sim, "_n_markers_"), negate = FALSE))
ids = as.numeric( which(ids) )
plots = bplapply(ids , function(id){
print(id)
p = get_plot(id , red_dim)
return(p)
}, BPPARAM = mcparam)
p = ggarrange(plotlist = plots , nrow = 12 , ncol = 12 , common.legend = T)
return(p)
}
anno = expand.grid(round_sim = c(1,2,3) , red_dim = c("all" , "hvgs"))
plots = lapply(1:nrow(anno) , function(i){
p = get_all_plots(as.numeric(as.character(anno$round_sim[i])) , red_dim = paste0("pca.", as.character(anno$red_dim[i])))
ggsave(filename = paste0(figures.dir, "umaps/" , as.character(anno$red_dim[i]), "/" , "all_" , as.character(anno$round_sim[i]), ".png"), plot = p, width = 17, height = 22)
})
sessionInfo()
## R version 4.2.3 (2023-03-15)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 22.04.2 LTS
##
## Matrix products: default
## BLAS: /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3
## LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.20.so
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
## [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] stats4 stats graphics grDevices utils datasets methods
## [8] base
##
## other attached packages:
## [1] tibble_3.2.1 dplyr_1.1.1
## [3] MetBrewer_0.2.0 viridis_0.6.2
## [5] viridisLite_0.4.1 ggpubr_0.6.0
## [7] ggplot2_3.4.2 stringr_1.5.0
## [9] miloDE_0.0.0.9000 miloR_1.6.0
## [11] edgeR_3.40.2 limma_3.54.2
## [13] scran_1.26.2 scuttle_1.8.4
## [15] splatter_1.22.1 geneBasisR_0.0.0.9000
## [17] BiocParallel_1.32.6 SingleCellExperiment_1.20.1
## [19] SummarizedExperiment_1.28.0 Biobase_2.58.0
## [21] GenomicRanges_1.50.2 GenomeInfoDb_1.34.9
## [23] IRanges_2.32.0 S4Vectors_0.36.2
## [25] BiocGenerics_0.44.0 MatrixGenerics_1.13.0
## [27] matrixStats_0.63.0 BiocStyle_2.26.0
## [29] rmarkdown_2.21
##
## loaded via a namespace (and not attached):
## [1] utf8_1.2.3 tidyselect_1.2.0
## [3] grid_4.2.3 combinat_0.0-8
## [5] munsell_0.5.0 ScaledMatrix_1.6.0
## [7] ragg_1.2.5 codetools_0.2-19
## [9] statmod_1.5.0 future_1.32.0
## [11] tester_0.1.7 withr_2.5.0
## [13] batchelor_1.14.1 colorspace_2.1-0
## [15] highr_0.10 knitr_1.42
## [17] ggsignif_0.6.4 pbmcapply_1.5.1
## [19] listenv_0.9.0 labeling_0.4.2
## [21] GenomeInfoDbData_1.2.9 mnormt_2.1.1
## [23] polyclip_1.10-4 optimParallel_1.0-2
## [25] farver_2.1.1 coda_0.19-4
## [27] parallelly_1.35.0 vctrs_0.6.2
## [29] generics_0.1.3 clusterGeneration_1.3.7
## [31] ipred_0.9-14 xfun_0.38
## [33] timechange_0.2.0 R6_2.5.1
## [35] doParallel_1.0.17 ggbeeswarm_0.7.1
## [37] graphlayouts_0.8.4 rsvd_1.0.5
## [39] locfit_1.5-9.7 pals_1.7
## [41] RcppZiggurat_0.1.6 bitops_1.0-7
## [43] cachem_1.0.7 DelayedArray_0.24.0
## [45] scales_1.2.1 nnet_7.3-18
## [47] ggraph_2.1.0 beeswarm_0.4.0
## [49] gtable_0.3.3 beachmat_2.14.2
## [51] globals_0.16.2 phangorn_2.11.1
## [53] tidygraph_1.2.3 timeDate_4022.108
## [55] rlang_1.1.0 systemfonts_1.0.4
## [57] scatterplot3d_0.3-43 splines_4.2.3
## [59] rstatix_0.7.2 yardstick_1.1.0
## [61] dichromat_2.0-0.1 broom_1.0.4
## [63] checkmate_2.1.0 reshape2_1.4.4
## [65] BiocManager_1.30.20 yaml_2.3.7
## [67] abind_1.4-5 backports_1.4.1
## [69] Rfast_2.0.7 lava_1.7.2.1
## [71] tools_4.2.3 bookdown_0.33
## [73] Augur_1.0.3 jquerylib_0.1.4
## [75] RColorBrewer_1.1-3 plyr_1.8.8
## [77] Rcpp_1.0.10 parsnip_1.1.0
## [79] sparseMatrixStats_1.11.1 zlibbioc_1.44.0
## [81] purrr_1.0.1 RCurl_1.98-1.12
## [83] rpart_4.1.19 cowplot_1.1.1
## [85] zoo_1.8-12 ggrepel_0.9.3
## [87] cluster_2.1.4 furrr_0.3.1
## [89] magrittr_2.0.3 magick_2.7.4
## [91] data.table_1.14.8 ResidualMatrix_1.8.0
## [93] lmtest_0.9-40 patchwork_1.1.2
## [95] evaluate_0.20 paleotree_3.4.5
## [97] gridExtra_2.3 scater_1.26.1
## [99] compiler_4.2.3 maps_3.4.1
## [101] htmltools_0.5.5 tidyr_1.3.0
## [103] expm_0.999-7 lubridate_1.9.2
## [105] DBI_1.1.3 tweenr_2.0.2
## [107] MASS_7.3-58.3 Matrix_1.5-4
## [109] car_3.1-2 wesanderson_0.3.6
## [111] cli_3.6.1 quadprog_1.5-8
## [113] gdata_2.18.0.1 parallel_4.2.3
## [115] metapod_1.6.0 gower_1.0.1
## [117] igraph_1.4.2 pkgconfig_2.0.3
## [119] RcppGreedySetCover_0.1.0 numDeriv_2016.8-1.1
## [121] recipes_1.0.5 foreach_1.5.2
## [123] vipor_0.4.5 bslib_0.4.2
## [125] ggcorrplot_0.1.4 hardhat_1.3.0
## [127] dqrng_0.3.0 XVector_0.38.0
## [129] prodlim_2023.03.31 digest_0.6.31
## [131] RcppAnnoy_0.0.20 phytools_1.5-1
## [133] fastmatch_1.1-3 uwot_0.1.14
## [135] DelayedMatrixStats_1.20.0 gtools_3.9.4
## [137] lifecycle_1.0.3 nlme_3.1-162
## [139] jsonlite_1.8.4 carData_3.0-5
## [141] BiocNeighbors_1.16.0 mapproj_1.2.11
## [143] fansi_1.0.4 pillar_1.9.0
## [145] lattice_0.21-8 survival_3.5-5
## [147] fastmap_1.1.1 plotrix_3.8-2
## [149] glue_1.6.2 png_0.1-8
## [151] iterators_1.0.14 bluster_1.8.0
## [153] ggforce_0.4.1 class_7.3-21
## [155] stringi_1.7.12 sass_0.4.5
## [157] textshaping_0.3.6 rsample_1.1.1
## [159] BiocSingular_1.14.0 future.apply_1.10.0
## [161] irlba_2.3.5.1 ape_5.7-1